********************************************************************************
* BACKGROUND
********************************************************************************
/*	This code analyzes how the COVID-19 pandemic affected removal outcomes for
	Chinese nationals in U.S. immigration courts and examines whether the 
	effect varied by immigration judge (IJ) characteristics. It begins by 
	plotting observed removal rates over time for Chinese and non-Chinese 
	nationals using margins and marginsplot. The code then runs a series of
	logistic regression models with triple interaction terms between Chinese
	nationality, the pandemic period, and IJ or court-level characteristics 
	such as race, gender, and court diversity. It calculates and compares the
	marginal effects of the pandemic across subgroups by using margins and 
	mlincom, focusing on differences between Chinese and non-Chinese respondents.
	For continuous court diversity measures, it evaluates the effect at the 
	median and at one standard deviation above the median. The script also 
	generates descriptive statistics for key variables. Finally, it constructs
	a new dataset from the model output and produces a figure that visually 
	compares the size and confidence intervals of pandemic effects across 
	judge subgroups*/

********************************************************************************
* READ IN DATA AND SET WORKING DIRECTORY TO SAVE OUTPUT FROM ANALYSIS
********************************************************************************

use pandemic_replication_data.dta, clear 
cd "output"
	
********************************************************************************
* FIGURE: OBSERVED REMOVAL RATES BY NATIONALITY PRE/POST PANDEMIC
********************************************************************************
	
	reg removal_decision i.p_comp_month##i.Chinese ///
		if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1
		margins p_comp_month#Chinese
		
*actual graph command 
	marginsplot, xdimension(p_comp_month) plotdimension(Chinese) ///
		recast(line) plot1opts(lcolor(gray) ///
		lwidth(medium) lpattern(dash)) ///
		plot2opts(lcolor(black) lwidth(medium)) ///
		recastci(rarea) ciopts(fcolor(%15) lcolor(%0)) ytitle(Removal Rate) ///
		yscale(lcolor(black) extend fextend) ///
		ylabel(, tlcolor(black)) xtitle("") ///
		xscale(lcolor(black) fextend) xline(722) ///
		xlabel(#5, format(%tmMonth_CCYY) tlcolor(black)) title("")
			
	graph export ///
		pandemic_CI_line.tif, ///
		width(3900) ///
		as(tif) replace	
	
********************************************************************************
* REGRESSION MODELS: PANDEMIC W CHINESE V OTHER COMPARISON
********************************************************************************

*categorical vars 
foreach var in ij_asian female hispanic_ij {
	
	logit removal_decision ///
		atty i.lang i.custody_n ///
		cube_rt_caseload i.IJ_FE_ID ///
		i.pandemic_MAR1##i.Chinese ///
		1.`var'#1.Chinese ///
		1.`var'#1.pandemic_MAR1 ///
		1.pandemic_MAR1#1.`var'#1.Chinese ///
		1.apt_dem_pres#1.Chinese ///
		1.apt_dem_pres#1.pandemic_MAR1 ///
		1.pandemic_MAR1#1.apt_dem_pres#1.Chinese ///
		if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1, cluster(nat_fe)

log using  pan_`var'.log, replace 
	
	margins `var'#Chinese, dydx(pandemic_MAR1) post
		mlincom, clear 
		mlincom (3-1), ///
		stat(est se p) add ///
		rowname("Non CN")		
		mlincom (4-2), ///
		stat(est se p) add ///
		rowname("CN")	
		mlincom (2-1), ///
		stat(est se p) add ///
		rowname(`var'==0)	
		mlincom (4-3), ///
		stat(est se p) add ///
		rowname(`var'==1)		
  		mlincom (4-2)-(3-1), ///
		stat(est se p) add ///
		rowname("DDD")

log close 
}

*effect of percent_asian and percent_hispanic among white judges 
foreach var in percent_asian_14 percent_hispanic_14 {
	
	logit removal_decision ///
		atty i.lang i.custody_n ///
		cube_rt_caseload i.IJ_FE_ID ///
		c.`var' ///
		i.pandemic_MAR1##i.Chinese ///
		c.`var'#1.Chinese ///
		c.`var'#1.pandemic_MAR1 ///
		1.pandemic_MAR1#c.`var'#1.Chinese ///
		1.apt_dem_pres#1.Chinese ///
		1.apt_dem_pres#1.pandemic_MAR1 ///
		1.pandemic_MAR1#1.apt_dem_pres#1.Chinese ///
		if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1 & ij_race=="White", cluster(nat_fe)

log using  pan_`var'.log, replace 
	
	sum `var' if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1 & ij_race=="White", d
		local median_val = `r(p50)'
		local SD = `r(sd)'
		local median_SD = `median_val' + `SD'
	
	margins Chinese, dydx(pandemic_MAR1)  ///
		at(`var'= (`median_val' `median_SD')) post 
		mlincom, clear 
		mlincom (3-1), ///
		stat(est se p) add ///
		rowname("Non CN")		
		mlincom (4-2), ///
		stat(est se p) add ///
		rowname("CN")			
		mlincom (2-1), ///
		stat(est se p) add ///
		rowname("Median")
		mlincom (4-3),  ///
		stat(est se p) add ///
		rowname("Median+SD")
  		mlincom (4-2) -(3-1), ///
		stat(est se p) add ///
		rowname("DDD")

log close 	
	}

*effect of percent female among male judges 
foreach var in  percent_female_14 {
	
	logit removal_decision ///
		atty i.lang i.custody_n ///
		cube_rt_caseload i.IJ_FE_ID ///
		c.`var' ///
		i.pandemic_MAR1##i.Chinese ///
		c.`var'#1.Chinese ///
		c.`var'#1.pandemic_MAR1 ///
		1.pandemic_MAR1#c.`var'#1.Chinese ///
		1.apt_dem_pres#1.Chinese ///
		1.apt_dem_pres#1.pandemic_MAR1 ///
		1.pandemic_MAR1#1.apt_dem_pres#1.Chinese ///
		if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1  & female==0, robust

log using  pan_`var'.log, replace 
	
	sum `var' if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1 & female==0, d
		local median_val = `r(p50)'
		local SD = `r(sd)'
		local median_SD = `median_val' + `SD'
	
	margins Chinese, dydx(pandemic_MAR1) ///
		at(`var'= (`median_val' `median_SD')) post 
		mlincom, clear 
		mlincom (3-1), ///
		stat(est se p) add ///
		rowname("Non CN")		
		mlincom (4-2), ///
		stat(est se p) add ///
		rowname("CN")	
		mlincom (2-1), ///
		stat(est se p) add ///
		rowname("Median")
		mlincom (4-3), ///
		stat(est se p) add ///
		rowname("Median+SD")
  		mlincom (4-2)-(3-1), ///
		stat(est se p) add ///
		rowname("DDD")
		
log close 	
	}

********************************************************************************
* DESCRIPTIVE STATISTICS
********************************************************************************

log using pandemic_descriptives.log, replace 

	sum removal_decision Chinese pandemic_MAR1 ///
		ij_asian hispanic_ij female ///
		percent_asian_14 percent_hispanic_14 percent_female_14 ///
		apt_dem_pres atty cube_rt_caseload ///
		English Spanish Other_Unknown_Lang ///
		Never_Detained Released Detained ///
		if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1
	
	mean removal_decision ///
		if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1 ///
		,over(Chinese pandemic_MAR1)

		log close 

********************************************************************************
* FIGURE: DIFFERENCES IN PANDEMIC EFFECTS   
********************************************************************************
*categorical vars 

foreach var in ij_asian hispanic_ijfemale {
	
	logit removal_decision ///
		atty i.lang i.custody_n ///
		cube_rt_caseload i.IJ_FE_ID ///
		i.pandemic_MAR1##i.Chinese ///
		1.`var'#1.Chinese ///
		1.`var'#1.pandemic_MAR1 ///
		1.pandemic_MAR1#1.`var'#1.Chinese ///
		if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1, cluster(nat_fe)
	
	margins `var'#Chinese, dydx(pandemic_MAR1) post
		mlincom, clear 
  		mlincom (4-2)-(3-1), ///
		stat(est ul ll ) add ///
		rowname("DD")
		gen `var'_e = _mlincom[1, 1] 
		gen `var'_ul = _mlincom[1, 2] 
		gen `var'_ll = _mlincom[1, 3] 
}

*effect of percent_asian and percent_hispanic among white judges 	
foreach var in  percent_asian_14 percent_hispanic_14 {
	
	logit removal_decision ///
		atty i.lang i.custody_n ///
		cube_rt_caseload i.IJ_FE_ID ///
		c.`var'##i.pandemic_MAR1##i.Chinese ///
		if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1 & ij_race=="White", cluster(nat_fe)

	sum `var' if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1 & ij_race=="White", d
		local median_val = `r(p50)'
		local SD = `r(sd)'
		local median_SD = `median_val' + `SD'
	
	margins Chinese, dydx(pandemic_MAR1) ///
		at(`var'= (`median_val' `median_SD')) post 
		mlincom, clear 
  		mlincom (4-2)-(3-1), ///
		stat(est ul ll) add ///
		rowname("DD")
		gen `var'_e = _mlincom[1, 1] 
		gen `var'_ul = _mlincom[1, 2] 
		gen `var'_ll = _mlincom[1, 3] 	
}	

*effect of percent female among male judges 
foreach var in percent_female_14 {
	
	logit removal_decision ///
		atty i.lang i.custody_n ///
		cube_rt_caseload i.IJ_FE_ID ///
		c.`var'##i.pandemic_MAR1##i.Chinese ///
		if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1 & female==0, robust

	sum `var' if include_m==1 &is_missing==0 & ///
		IJ_pre_post_pan==1 & female==0, d
		local median_val = `r(p50)'
		local SD = `r(sd)'
		local median_SD = `median_val' + `SD'
		
	margins Chinese, dydx(pandemic_MAR1) ///
		at(`var'= (`median_val' `median_SD')) post 
		mlincom, clear 
  		mlincom (4-2)-(3-1), ///
		stat(est ul ll) add ///
		rowname("DD")
		gen `var'_e = _mlincom[1, 1] 
		gen `var'_ul = _mlincom[1, 2] 
		gen `var'_ll = _mlincom[1, 3] 
}	

*store output from margins 
	frame put /// 
	ij_asian_ll hispanic_ij_ll female_ll ///
	percent_asian_14_ll percent_hispanic_14_ll percent_female_14_ll ///
	ij_asian_ul hispanic_ij_ul female_ul ///
	percent_asian_14_ul percent_hispanic_14_ul percent_female_14_ul ///
	ij_asian_e hispanic_ij_e  female_e ///
	percent_asian_14_e percent_hispanic_14_e percent_female_14_e, ///
	into(fig_pan)
 
*make frame/labels  to plot stats 
	frame change fig_pan
	duplicates drop 
	xpose, clear varname 
	rename (v1) (Pan_effect)
	generate str var_id = substr( _varname , 1, strlen( _varname) - 2)
	gen yes_no=substr(_varname,-2,.) 
	
	replace var_id ="6. Percent Female Court" if strpos( var_id, "percent_female")
	replace var_id ="5. Percent Hisp Court" if strpos( var_id, "percent_hispan")
	replace var_id ="3. Female" if strpos(var_id, "female")
	replace var_id ="1. Asian" if strpos( var_id, "ij_asian")
	replace var_id ="2. Hispanic" if strpos( var_id, "hispanic_ij")
	replace var_id ="4. Percent Asian Court" if strpos( var_id, "percent_asian")
	
	gen dif = Pan_effect if yes_no=="_e"
	gen ul = Pan_effect if yes_no=="ul"
	gen ll = Pan_effect if yes_no=="ll"
   
	bysort var_id : egen Difference = max(dif)
	bysort var_id : egen U_CI = max(ul)   
	bysort var_id : egen L_CI = max(ll)   
	keep var_id  Difference U_CI L_CI
	duplicates drop 
   
	encode var_id, gen(var_idn)
	label define graph_ab_difn ///
		1 "Asian IJ" ///
		2 "Hispanic IJ" ///
		3 "Female IJ" ///
		4  "Percent Asian Court +SD" ///
		5  "Percent Hispanic Court +SD" ///
		6 "Percent Female Court +SD" 
	
	label values var_idn graph_ab_difn

*actual plot 
	twoway (rcap U_CI L_CI var_idn, sort) ///
	(scatter Difference var_idn, mcolor(black) ///
	msize(medium) msymbol(circle) mlabel(var_idn)  mlabcolor(black) ///
	mlabposition(9) mlabangle(vertical) mlabgap(small)), ///
	ytitle("Differences in Pandemic Effects" "for Chinese vs. Non-Chinese Nationals") ///
	yscale(noline) yline(0) ///
	ylabel(, ticks tlcolor(black)) xtitle("") ///
	xlabel(.85 1 2 3 4 4.15, nolabels noticks nogmin nogmax) ///
	legend(off) plotregion(lcolor(black) ilcolor(none))
	
	graph export ///
	pandemic_dif.tif, ///
	width(3900) replace 